Single-Cell Dissection of the Immune Response After Acute Myocardial Infarction

BACKGROUND: The immune system’s role in ST-segment–elevated myocardial infarction (STEMI) remains poorly characterized but is an important driver of recurrent cardiovascular events. While anti-inflammatory drugs show promise in reducing recurrence risk, their broad immune system impairment may induce severe side effects. To overcome these challenges, a nuanced understanding of the immune response to STEMI is needed. METHODS: For this, we compared peripheral blood mononuclear single-cell RNA-sequencing (scRNA-seq) and plasma protein expression over time (hospital admission, 24 hours, and 6–8 weeks post-STEMI) in 38 patients and 38 controls (95 995 diseased and 33 878 control peripheral blood mononuclear cells). RESULTS: Compared with controls, classical monocytes were increased and CD56dim natural killer cells were decreased in patients with STEMI at admission and persisted until 24 hours post-STEMI. The largest gene expression changes were observed in monocytes, associating with changes in toll-like receptor, interferon, and interleukin signaling activity. Finally, a targeted cardiovascular biomarker panel revealed expression changes in 33/92 plasma proteins post-STEMI. Interestingly, interleukin-6R, MMP9 (matrix metalloproteinase-9), and LDLR (low-density lipoprotein receptor) were affected by coronary artery disease–associated genetic risk variation, disease status, and time post-STEMI, indicating the importance of considering these aspects when defining potential future therapies. CONCLUSIONS: Our analyses revealed the immunologic pathways disturbed by STEMI, specifying affected cell types and disease stages. Additionally, we provide insights into patients expected to benefit most from anti-inflammatory treatments by identifying the genetic variants and disease stage at which these variants affect the outcome of these (drug-targeted) pathways. These findings advance our knowledge of the immune response post-STEMI and provide guidance for future therapeutic studies.


Introduction
Acute myocardial infarction (MI) is one of the leading causes of mortality and morbidity globally. The immune system plays an important role in the pathophysiology of the disease from the initial formation and rupture of atherosclerotic plaques, all the way up to inflammation in response to the MI (1-3).
Systemically studying and targeting the inflammatory response has remained challenging due to its complexity and severe side effects, as was seen in the CANTOS trial (interleukin-1ß, IL-1ß) and ASSAIL-MI (interleukin 6 receptor, IL6R) trial (4,5). Identifying specific immune cell subpopulations and underlying processes that are specifically disturbed after MI could open the way for more targeted treatment that is expected to cause fewer side effects (6,7).
However, until recently, analyzing the characteristics of specific cell types was often done at bulk level and using predefined marker genes, not providing an unbiased, comprehensive overview of the immune system (8). Recent developments in the field of single-cell RNA sequencing (scRNA-seq) now enable us to unbiasedly study 100,000s of individual cells simultaneously on a transcriptome-wide level (9)(10)(11). This technology enables us to gain new insights into the inflammatory response after a MI by flexibly and concurrently mapping MI-induced changes in cell type composition, gene expression level and downstream pathways at various cellular resolutions.
To study MI at the single-cell level, it is essential to first understand the normal physiological variation among the general population. Such groundwork has been generated both for the heart (12) and the circulating immune cells (10). Various aspects in the pathophysiology of myocardial infarction have now been dissected at the single-cell level, from atherosclerotic aortas precluding the disease (13,14) to cardiac neovascularization by resident endothelial cells after the infarction (15). More recently, single-cell multi-omic maps (combining spatial transcriptomics, single-cell expression and chromatin accessibility profiling) have been generated from human heart tissue after MI (16). This  Figure S1). Moreover, we show that several important plasma proteins during the disease are affected by a combination of genetics, disease status or phase of the disease. Altogether, the current study increases the resolution by which the immune system during and after MI has been dissected and emphasizes the importance of taking person-and disease-related characteristics into account to fully grasp the observed underlying molecular changes.
All rights reserved. No reuse allowed without permission.
(which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.
The copyright holder for this preprint this version posted May 2, 2023.

Study population
All patients presenting with a first STEMI upon admission at the Heart Catherization center at the University Medical Center Groningen (UMCG) were enrolled in the study between January 1st 2018 and November 30th 2019. Inclusion criteria were adults (>18 years) that had a STEMI, a primary percutaneous intervention (PCI) with implantation of at least one stent with a diameter of at least 3 mm resulting in thrombosis in myocardial infarction (TIMI) flow grade 2 or 3 post PCI, and which showed less than 6 hours symptoms before undergoing PCI. Major exclusion criteria were previous MI, medical history of diabetes, inflammatory disease or malignancies, medication affecting inflammation and clemastine or desloratadine use during intervention. Using these criteria, we initially enrolled 88 STEMI patients in this study, of which 42 were later excluded due to missing follow-up time points (t24h, t8w).
Further 8 patients were excluded in the follow-up due to presence of exclusion criteria that were not yet known upon admission (Figure 1). The study was part of CardioLines, a single-center observational biobank aimed to study potential factors related to success or failure of diagnosis and treatment both from a patient as well as a medical perspective and was approved by the ethics committee of the UMCG, document number METC UMCG 2012/296 (17). Informed consent was obtained of all patients. scRNAseq data from 38 age-and sex-balanced participants from the LifeLines DEEP cohort were included as a control group (18,19).

Clinical parameters
Upon inclusion of STEMI patients whole blood samples were collected at three different timepoints: during admission at the heart catheterisation center (t0), 24 hours (t24h, acute phase) and 6-8 weeks (t8w, chronic phase) later (Figure 1). Standard laboratory assessment of the blood was performed upon admission and routine physical parameters were recorded ( Table 2). Titres of creatine kinase (CK), 10 myocardial band of CK (CK-MB) and troponin T were routinely measured at 3h, 6h, 9h, 12h, 24h and 48h after admission.

Isolation and preparation of PBMCs
PBMCs were isolated and stored as previously reported (9). In short, for each donor and timepoint, whole blood was collected in two 10 ml EDTA-vacutainer (BD) tubes. Within 2 hours after collection, PBMCs were isolated from the blood using cell preparation tubes with sodium heparin (BD). PBMCs were kept in RPMI1640 supplemented with 50 µg/mL gentamicin, 2 mM l-glutamine, and 1 mM pyruvate. Isolated PBMCs were cryopreserved in RPMI1640 containing 40% FCS and 10% DMSO. Within one year, PBMCs were further processed for scRNA-seq analysis. Cells were thawed in a 37°C water bath until almost completely thawed, after which the cells were washed with pre-warmed (37°C) RPMI1640.
After washing, cells were resuspended in pre-warmed RPMI1640 medium and incubated for 1 h in a 45° slant rack at 37°C in a 5% CO2 incubator. After resting, cells were washed twice in medium supplemented with 0.04% bovine serum albumin. Cells were counted using a hemocytometer and cell viability was assessed by Trypan Blue assay.

Single-cell library preparation and sequencing
Fifteen sex-and timepoint-balanced sample pools were prepared aiming for 1,000 cells per patient for a total of 8 patients per pool. Single cells were captured with the 10X Chromium controller (10x Genomics) according to the manufacturer's instructions (document CG00026) and as described earlier (20 (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

Alignment and initial processing of sequencing data
CellRanger v3.0.2 software with default settings was used to demultiplex the sequencing data, generate FASTQ files, align the sequencing reads to the hg19 reference genome, filter cell and UMI (unique molecular identifier) barcodes, and count gene expression per cell. Genotypes of the samples were generated using the Infinium Global Screen Array-24 kit (v1 for the samples processed with 10X v2 chemistry and v3 for the samples processed with 10X v3 chemistry) and were phased using Eagle v2.4 and imputed with the HRC 1000G Phase-3 v5 (hg19) reference panel using the Michigan Imputation Server (Minimac v3) (21). Genotype information of the controls was processed as described in the Lifelines DEEP cohort paper (18). Doublet detection and sample assignment were done using Souporcell (22), confining ourselves to only use exonic single nucleotide polymorphisms (SNPs) with a MAF≥0.05.
The detected doublets were discarded. Then, samples were assigned by correlating the genotypes of the clusters in each 10X lane with the genotype information of the samples in the respective lane, and assigning to the sample with the highest correlation.

Data preparation and quality control
Version 4.0 of the R package Seurat was used for data preprocessing (23)(24)(25). We performed dimensional reduction and clustering to integrate four datasets (v2 and v3 10X chemistry for STEMI patients and controls). For this, the Pearson residuals were used resulting from running SCTransform in each dataset using Seurat's SCTintegration method (26). The first 30 principal components were used for cell clustering using Seurat's FindNeighbors (dims = 1:30) and FindCluster function (default parameters, All rights reserved. No reuse allowed without permission. (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.
The copyright holder for this preprint this version posted May 2, 2023.
12 resolution 1.2) and an UMAP plot was used to visualize this. All genes that were not detected in at least 3 cells were removed. Cells with a high percentage of reads mapping to mitochondrially-encoded genes (>8% for 10X V2 and >15% for 10X V3) were discarded, as this can be a marker of poor-quality cells (27). Also, cells expressing ≤200 genes were considered outliers and discarded. Furthermore, cells expressing ≥ 10 UMIs of the HBB gene were discarded, indicative of erythrocytes. Finally, all cells that were marked as doublets or inconclusive by the SoupOrCell method were discarded (22). In total, 129,873 cells were used for downstream analysis consisting of 95,995 cells from STEMI patients and 33,878 cells from ageand sex-balanced controls ( Table S1).

Cell type abundance
Cell type proportion differences across STEMI timepoints and differences with the controls were determined using the permutation-based method described in Farbehi et al. (29). The number of cells of All rights reserved. No reuse allowed without permission.
(which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.
Pairwise comparisons were made across all combinations of STEMI timepoints and the controls using a W-value of 0.05. Resulting p-values were Holm-Bonferroni corrected, with a corrected p-value < 0.05 being considered statistically significant.

Differential gene expression and pathway analysis
For differential gene expression (DE) analysis library-size normalization was performed by scaling the total UMI gene expression of each cell barcode to 10,000, after which the gene expression data was logtransformed (LogNormalize in Seurat). For each cell type, pairwise DE analysis was performed using MAST comparing everything to t0 (controls, 24 hours and 6-8 weeks) (30). Due to strong technical differences in gene expression profiles between the samples that were prepared with the 10X V2 or 10X V3 protocol, DE analysis was performed separately for the samples belonging to either version of the protocol. A meta-analysis, considering only genes tested in both datasets, was performed across the V2 and V3 results using the Fisher's Combined Probability Test. To correct for multiple testing, Bonferroni correction was applied and adjusted p-values <0.05 were considered statistically significant. The up-and downregulated DE genes were used for separate pathway analyses using the ToppFun functional enrichment tool, selecting for the REACTOME database to identify relevant pathways (31).

Cell-to-cell communication analysis
We predicted potential cell-to-cell interactions using NicheNet analysis v1.1.0 (32). Using a background set of DE genes, NicheNet predicts potential ligand-receptor interactions and predicts whether specific cell-cell interactions led to differences in downstream gene-expression. We modified NicheNet's FindMarkers wrapper to use the MAST implementation in Seurat in order to determine DE ligands and receptors, as well as genes downstream of cells potentially interacting with the ligands. Differential 14 ligand and downstream expression changes were confined to genes that were present in at least 10% of the cells. The absolute log fold change between the conditions needed to be at least 0.1 to be considered for testing. We assessed potential ligand-receptor interactions between pairs found in the Omnipath database (33). The Nichenet supplied ligand-receptor network data was used to determine what could be considered to be cell communication patterns. Predicted ligand-receptor interactions were filtered for ligands that were shown to have an effect on downstream gene expression. Ligandreceptor connections between cell types were counted and visualized in circle plots.

Protein Olink analysis
Plasma from the 38 STEMI patients was collected at each of the three timepoints (t0, t24h, t8w), but missing one sample at t8w (n=37), and stored at -80°C. Plasma samples were analyzed for 92 cardiovascular disease related protein biomarkers by Olink Proteomics (Olink Target 96 Cardiovascular III panel, Uppsala Sweden). The data received was QCed by excluding all proteins with a missing data frequency of >15% and samples with an Olink warning QC outcome from analysis. All samples with a low detection rate were given the limit of protein detection value for the specific protein. Using these QC parameters, none of the proteins had to be excluded. Protein quantitative trait loci (pQTL) analysis was done using our in-house eQTL mapping pipeline [https://github.com/molgenis/systemsgenetics/wiki/QTL-mapping-pipeline]. pQTLs were mapped per time point, only testing the SNPs within a 100k window of the transcription start site of the gene transcribing the tested protein.

Relating clinical variables to cell type abundance, gene and protein expression
Cell type proportions and gene and protein expression profiles were associated with peak CK-MB values in STEMI patients. For gene expression associations we used MAST, a two-part generalized linear model All rights reserved. No reuse allowed without permission.
(which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.
The copyright holder for this preprint this version posted May 2, 2023. ; https://doi.org/10.1101/2023.05.02.23289370 doi: medRxiv preprint 15 on the log-transformed count matrix, and correcting for age, sex, 10X chemistry and 10X lane. For cell type proportion and protein expression associations, we first conducted a centered-log-ratio or logtransformation, respectively, after which linear regression was conducted, correcting for sex and age.

pQTL analysis
For STEMI patients individual Olink plasma protein levels were linked to specific genetic variants (pQTL mapping) at each of the 3 timepoints (t0, t24h, t8w) using the eQTL mapping pipeline (https://github.com/molgenis/systemsgenetics). We limited ourselves to the 241 variants that were previously found to show a genome-wide significant GWAS signal for coronary artery disease (CAD) (34). The significant pQTLs in STEMI patients in at least one timepoint were then assessed in controls by mapping these pQTLs in 1142 individuals from the general population cohort Lifelines Deep, in which the same Cardiovascular Olink panel was previously measured (35) and genotype data was available.
Genetic variants were converted to Trityper format using GenotypeHarmonizer (https://github.com/molgenis/systemsgenetics). To assess whether the pQTL effect strength was influenced by the timepoint after the STEMI (t0, t24h, t8w), we then conducted an interaction analysis by fitting a generalized linear model correcting for age and sex as fixed effects, and the participant as a random effect.
All rights reserved. No reuse allowed without permission.
(which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

Statistical analysis
Discrete variables were represented as frequencies and percentages. Continuous variables with a normal distribution were summarized as mean ± standard deviation and if skewed were represented as median with interquartile range. Spearman correlation was used for non-normally distributed variables.
For DE gene and protein analysis, and their correlation to peak CK-MB values, multiple testing correction was performed using Holm-Bonferroni. We considered a corrected P-value <0.05 statistically significant.
Statistical analysis and the figures were created using R Core Team (2020).

Population characteristics
To dissect the immune response during and after a MI, PBMCs and plasma from 38 patients with a first STEMI were collected at the time of hospital admission (t0), and 24h (t24h) and 6-8 weeks (t8w) after the event (Figure 1). The patient characteristics during hospital admission are described in Table 1. On average, STEMI patients were 60±11 years, had a median BMI of 26.4 (SD 3.5) kg/m 2 and the majority were males (84%). Most patients presented with one stenosed vessel (55%), and in 39% of the patients the affected vessel was completely occluded (TIMI 0 flow). These 38 STEMI patients were compared to 38 age-and sex-balanced controls from the LifeLines DEEP cohort whose scRNA-seq data was previously published (19). These controls showed comparable age-(59.2 years) and sex-characteristics (84% male), but a slightly lower BMI (25.4 (SD 3.3) kg/m 2 ) than the STEMI patients ( Table 1).
All rights reserved. No reuse allowed without permission.
(which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.
The copyright holder for this preprint this version posted May 2, 2023. ; Figure 1 (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

Single-cell profiling of immune cells in STEMI patients and controls
The collected PBMCs from the STEMI patients at three different timepoints post-STEMI (t0, t24h, t8w) and controls from the LifeLines cohort (19), were then used for 10X Genomics scRNA-seq analysis. A combination of 10X Genomics v2 and v3 chemistry reagents were used to capture an average of 842 cells per individual per condition (v2: 831 genes/cell, v3: 1,548 genes/cell) in STEMI patients, and an average of 891 cells per individual in controls (v2: 1,012 genes/cell, v3 1,931 genes/cell) ( Table S2).
Souporcell was used to demultiplex and to identify the doublets coming from different individuals (22).
Souporcell clusters were correlated to individuals by correlating the cluster genotypes with donor genotypes using Spearman correlation. This revealed that by mistake we omitted loading cells for one donor at the t8w condition and that on average, 10.0% of the cells were assigned as doublet or inconclusive. Due to technical differences between v2 and v3 chemistry, determination of quality control (QC) thresholds and analyses were performed separately per chemistry. Low quality cells were excluded, leaving 129,873 cells (95,995 diseased and 33,878 control) in the final dataset used for analyses (see Methods). UMAP dimensionality reduction and KNN-clustering was then applied on the normalized, integrated count data, allowing the identification of cell types.

Cell type abundance
We identified 10 major cell types in the PBMC fraction of STEMI patients and controls, including the B, CD4 + T, CD8 + T, DC, HSPC, monocyte, NK, plasmablast, platelet and other T cells (Figure 2a). At higher resolution, these major cell types could be further splitted in 29 minor cell populations, including different subsets of T cells (regulatory T cells (Treg), CD4 + and CD8 + T-naive, CD4 + and CD8 + T central memory (TCM), CD4 + and CD8 + T effector memory (TEM), CD4 + cytotoxic (CTL), CD4 + and CD8 + All rights reserved. No reuse allowed without permission.
(which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.
The copyright holder for this preprint this version posted May 2, 2023. To conduct a robust cell type composition analysis between the different conditions (STEMI patients at t0 compared to controls, STEMI patients over time), we initially focused our analyses on the 6 most abundant cell types (≥1% of the total PBMCs in each of the tested conditions): B cells, CD4 + T cells, CD8 + T cells, DCs, monocytes and NK cells (Figure 2c, S1, Table S2, S3). After which, we zoomed in at the higher resolution subcell type compositional changes to define the subtypes that contributed to these changes in composition (Figure S3).
(which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

Gene expression profiles of circulating immune cells in STEMI patients and controls
To assess gene expression differences in STEMI patients at time of hospital admission versus controls, we performed DE gene analysis on the six major cell types using MAST (Table S4). Of all cell types, All rights reserved. No reuse allowed without permission.
(which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.  (Figure 3a). The majority of the DE genes were uniquely identified in only one major cell type (Figure 3b). Looking at the monocytes, the cell type with most total and unique DE genes, we observed pathway enrichment among upregulated genes in STEMI patients involving various pro-inflammatory pathways, including IFN signaling, interleukins (IL-1, IL-6) and all toll-like receptors ( Table S5).
Subsequent comparison of the gene expression levels of STEMI patients at hospital admission compared to later time points revealed the largest amount of DE genes in NK cells and in monocytes, respectively 24h and 6-8 weeks after STEMI (Figure 3c, Table S4). Specifically looking at the pathways underlying the upregulated genes at these cell type-timepoint combinations, we observed that many of the NK DE genes at t24h were involved in interactions with platelets (activation, degranulation, response to platelet cytosolic Ca2+), whereas many of the monocyte DE genes at t8w were involved in proinflammatory signaling (IFN-gamma, IL-1) (Table S5). Remarkably, all major types showed a decrease in DE genes going from 24h to 6-8 weeks after STEMI, except for the monocytes (Figure 3c). And as before in comparison to controls (Figure 3b), monocytes also showed the largest proportion of unique DE genes over the course of the disease (Figure 3d). Together with, but independently of, the observed cell type composition changes (Figure 2c), these results indicate that both the monocytes and NK cells seem to play a pivotal role in the immunological changes that are occurring upon a STEMI. All rights reserved. No reuse allowed without permission.
(which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

Figure 3. Differential expressed genes and pathways in STEMI patients and in controls.
25 to-cell communication tool that identifies potential cell-cell interactions by assuming that ligandreceptor interactions can be predicted based on the expression of the ligand in one cell type and downstream gene expression changes of a known ligand-receptor interaction in another cell type (Table   S6). This analysis revealed that in the acute phase of the disease (t0-t24h) in-and outcoming communication was balanced among the cell types (Figure 4a). In contrast, in the chronic phase of the disease (t24h-t8w) the T cells and DCs were mostly sending communication, while monocytes and B cells were mostly receiving it (Figure 4b). Although from the comparable number of incoming interactions for a specific cell type it may seem that the involved ligands are mostly shared among cell types, most of these ligands were actually uniquely involved in the communication by just one cell type (Figure 4c, 4d).
Looking at the number of ligand-receptor interactions that were associated with changes in downstream signaling, we observed both NK and CD8 + T cells to be the largest communicators during the acute phase (Figure 4a), whereas monocytes were the main communicators during the chronic phase of the disease (Figure 4b). These results align well with our DE analyses, which also indicated the largest impact on expression changes of STEMI on NK cells early and on monocytes later in the disease.
All rights reserved. No reuse allowed without permission.
(which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

Protein expression profiles in STEMI patients over time
All rights reserved. No reuse allowed without permission.
(which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.
The copyright holder for this preprint this version posted May 2, 2023. ; https://doi.org/10.1101/2023.05.02.23289370 doi: medRxiv preprint While our previous analyses at the mRNA level help us to disentangle the molecular processes that immune cells undergo after STEMI, changes in plasma protein levels may provide us additional insights on the systemic consequences of a STEMI and can be more easily monitored in clinical follow-up of patients. Therefore, next, using the Olink Target 96 Cardiovascular III panel, we analyzed the plasma levels of 92 proteins that are known or exploratory human cardiovascular and inflammatory markers in each of the three timepoints (t0, t24h, t8w) within STEMI patients.  Table S7). Both upregulation of NT-proBNP and the soluble form of IL1RL1, a member of the IL-1 receptor family, are well known proteins that are upregulated in response due to increased wall stress during STEMI (36). Moreover, both are independent predictors of heart failure and cardiovascular death (37,38). CHI3L1, also known as YKL-40, is an extracellular matrix protein which is involved in several stages of CAD. It induces maturation of monocytes to macrophages in the early phase of atherosclerosis, but also promotes atherosclerotic progression and complications like plaque rupture (39). So, the upregulation of this protein might be due to the cause of MI (i.e. plaque rupture), whereas NT-proBNP and IL1RL1 reflect the mechanical stress to the cardiomyocytes due to MI.
(which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.
The copyright holder for this preprint this version posted May 2, 2023. ; https://doi.org/10.1101/2023.05.02.23289370 doi: medRxiv preprint 28 acute and chronic phase may be an indication of effective treatment by PCI after STEMI, resulting in restoration of blood flow and eventually preservation of heart function.
The top three downregulated proteins during the chronic phase were SPON1 (mean LFC=-1.2, p=1.3x10 -9 ), MPO (mean LFC=-0.36, p=2.0x10 -6 ) and AZU1 (mean LFC=-0.35, p=5.6x10 -4 ). Similar to the acute phase, both SPON1 and AZU1 continued to be downregulated. Additionally, MPO was found downregulated, a protein involved in the acute innate inflammatory response which contributes to plaque destabilization through local oxidative tissue injury (47,48). Together, these changes seem to indicate the long-term restoration of heart function and blood flow by PCI.
All rights reserved. No reuse allowed without permission.
(which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.
The copyright holder for this preprint this version posted May 2, 2023.

Relating Peak CK-MB as proxy for biochemical infarct size
All rights reserved. No reuse allowed without permission.
(which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.
The copyright holder for this preprint this version posted May 2, 2023. ; https://doi.org/10.1101/2023.05.02.23289370 doi: medRxiv preprint 30 STEMI patients can already widely differ from each other when they are presented with symptoms in the clinic. While the patient inclusion and exclusion criteria of our study (Methods) reduced the amount of clinical variation in our study (Table 1), many aspects remain that could underlie the observed donorto-donor variation in the measured molecular phenotypes, including age, sex, infarct size and genetics.
To assess whether variation in infarct size may correlate with the individual molecular data layers, we used plasma peak CK-MB levels as proxy. In clinical practice, plasma CK-MB levels are measured regularly after STEMI to monitor myocardial damage. This biomarker reaches a peak within 24 hours, and the peak value is an estimation of the infarct size and predicts left ventricle dysfunction (49). First, we used plasma peak CK-MB levels to assess whether PBMC cell type proportions and monocyte gene expression levels were associated with biochemical infarct size. We focussed on gene expression of the monocytes as this was the cell type with the largest amount of total and unique DE genes. However, neither cell type proportions, nor monocyte gene expression levels showed a significant association with peak CK-MB at any of the tested timepoints (t0, t24h or t8w). Next, we assessed whether peak CK-MB values were associated with inflammatory plasma protein levels in the acute and chronic phase after STEMI using simple linear regression. Of all 92 proteins, only NT-proBNP was positively associated with peak CK-MB in the acute phase (t24h: r 2 adjusted=0.49, p=2.8x10 -3 ) ( Table S8). NT-proBNP, a marker for myocardial wall stress, is commonly measured in clinical practice to assess left ventricular dysfunction (50,51). As there is a direct relationship between the level of myocardial damage and the amount of left ventricular dysfunction, NT-proBNP and CK-MB levels are expected to be correlated during the time CK-MB peaks (~24h after STEMI) (36,52).

QTL analysis
Beyond infarct size, genetic variation among patients, can be an important contributor to the observed differences in molecular response among patients. Therefore, it is important to take genetic variation All rights reserved. No reuse allowed without permission.
(which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.
The copyright holder for this preprint this version posted May 2, 2023. ; https://doi.org/10.1101/2023.05.02.23289370 doi: medRxiv preprint 31 into account when assessing a patient's disease course and defining patient-tailored treatment. As the measured plasma protein levels provide a pool of potential direct therapeutic targets, we assessed the effect of genetic variation on these (Table S9). To provide a direct clinical link, we specifically focussed on those genetic variants that were previously found to be associated with CAD risk (34). For three out of 92 plasma proteins, we detected a significant pQTL in the STEMI patients in at least one of the three timepoints (t0, t24h, t8w). For these three pQTLs, we assessed their disease-specificity by comparing pQTL effect sizes with those in a control cohort of the general population. This revealed one pQTL effect that was significant in all conditions and whose effect size was comparable among both patients (independent of time point) and controls: SNP rs6686750 affecting IL6R plasma protein levels (STEMI patients: genotype β=0.42, p=1.6x10 -5 , Figure 6a). The two other pQTLs were disease-or even condition-specific: SNP rs8124182 affecting MMP9 plasma protein levels was found specifically in STEMI patients (STEMI patients: genotype β=-0.68, p=0.016, Figure 6b) and SNP rs10422256 affecting LDLR plasma protein levels was found only in STEMI patients at t24h (STEMI patients: genotype x t24h interaction β=-0.29, p=0.022, Figure 6c). These results indicate the importance of taking genotype into consideration, when studying the relationship between STEMI and molecular phenotype.
All rights reserved. No reuse allowed without permission.
(which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

Discussion and Conclusion
In this study, we have provided an unbiased longitudinal overview of the single-cell immune response in STEMI patients and compared this to age-and sex-balanced controls. We first mapped changes in the cell type composition and gene and plasma protein expression level. We then determined the overarching pathways that were enriched and the cell types being involved in these changes as a result of cell-to-cell communication. Thereafter, we determined whether these changes were correlated with an individual's (biochemically-measured) infarct size, genotype or disease phase. Together, all these analyses indicated an important role for the classical monocytes and NKdim cells during STEMI, both in terms of cell type compositional and gene expression changes.
Our observation that monocytes play a pivotal role in STEMI is in line with a previous scRNA-seq study of 10 STEMI patients (with and without plaque rupture) that were analyzed at time of hospital admission (53). While this study did not observe any differences in cell type composition between STEMI patients with or without plaque rupture, they did observe that most cell-cell communication occurred from and to monocytes. Our study now adds additional insights, because of the longitudinal character, the added layer of protein information, the increased sample size, and the availability of controls to compare with. The recently published spatial multi-omic map of human MI, in which cardiac structure and specific cardiac cells were studied in spatial context, underlined the role of immune cells in cardiac repair following STEMI (16). While that study used cardiac tissue, and we studied the response of circulating immune cells, it remains clear that gaining more insights into the immune response after STEMI is key for future drug therapies. All rights reserved. No reuse allowed without permission.
(which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.
The copyright holder for this preprint this version posted May 2, 2023. ; https://doi.org/10.1101/2023.05.02.23289370 doi: medRxiv preprint One such insight relevant for future drug therapies that our study has now provided, is in relation to the use of IL-1ẞ blockers in MI. Our study showed chronic upregulation of the IL-1 signaling pathway specifically in the monocytes of STEMI patients. The CANTOS trial showed that IL-1ẞ inhibition with canakinumab led to lower cardiovascular event rates, but also higher fatal infections (4). Our insights into the IL-1 pathway could aid in targeting specifically monocytes, and thereby potentially preventing the occurrence of serious side effects. Furthermore, a study in IL-1R knock-out mice showed that these mice, which cannot respond to IL-1, were protected from hypertrophy, fibrosis, and diastolic dysfunction (56). This may indicate that the observed chronic upregulation of the IL-1 signaling pathway after STEMI reflects an already ongoing process of developing heart failure and eludes the potential to broaden the use of IL-1 inhibition for both atherosclerotic events and fibrotic processes.
Additionally, our study indicated the importance of considering both genetic variation and disease phase when assessing the molecular consequences of a STEMI, as several of the protein expression levels were affected by a combination of these parameters. Importantly, by combining molecular layers, our integrative analysis revealed insights in several therapeutic targets for which drugs are already available and are currently tested in clinical trials.
The first example for a promising drug target in the context of CVD, involves targeting IL-6 signaling. IL-6 signaling can be activated through three different modes of signaling, each having their own consequences: classical signaling, trans-signaling and trans-presentation (54). The antiinflammatory classical signaling is activated by membrane-bound IL6R and gp130 (also known as IL6ST) within the same cell, the pro-inflammatory trans-signaling is regulated by soluble IL-6R and membranebound gp130, and the Th17-cell promoting trans-presentation signaling is induced by membrane-bound IL-6R (in DCs) and gp130 (in Th17 cells). Through cell-cell communication analyses with NicheNet (Table   S6), we found that both the classical and trans-presentation pathways were differentially active during the first 24h of STEMI. While our pQTL analysis indicated that the trans-signaling pathway was affected All rights reserved. No reuse allowed without permission.
(which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.
Each of these IL-6 pathways can be targeted using different drugs (54), and based on our findings it seems that these drugs may be most beneficial during specific stages of the disease (first 24h after STEMI: for classical and trans-presentation IL-6 signaling, using drugs like anti-IL-6 (54,55) or anti-IL-6R) or in donors with a specific genetic background (patients with the GG genotype at SNP rs6689206 have the highest expression of IL-6R and thereby will naturally activate IL-6 trans signaling the most, and as such are expected to benefit most of treatments targeting specifically this pathway, e.g. Olamkicept, soluble gp130) (Fig 6A). Another interesting drug target for which we now provide additional insights is the LDLR. LDLR plasma levels act as scavengers for plasma LDL and can be indirectly targeted by PCSK9 inhibitors, which inhibit the degradation of membrane LDLR (44). As such, they are widely used to decrease LDL cholesterol and reduce the risk of cardiovascular events (44). Our DE protein analyses revealed that PCSK9 levels increased after STEMI, while our pQTL analysis revealed that 24h after STEMI individuals with the GG genotype (at SNP rs10422256) have a lower plasma protein level of LDLR independent of their PCSK9 levels and may therefore have higher free plasma LDL leading to increased cardiovascular disease risk (Fig. 6C). Together, this indicates that individuals with the GG genotype at SNP rs10422256 may benefit the most from PCSK9 inhibitor treatment, especially in the first 24 hours of the disease. We expect that implementing such genetic and phase of the disease information will be important to maximize benefit of such patient-tailored therapies in future clinical practice.
Unfortunately, when focusing on how each of the molecular layers could be contributing to or being the consequence of biochemically-measured infarct size, we could not identify any direct relationship except for plasma NT-proBNP levels -a known marker for heart failure (50,51). We suspect this lack of association might have been the result of assessing each molecular parameter separately or not considering all potentially relevant clinical and donor variables (age, sex, smoking status, genetics, etc.) into the model 57 . To enable such a model considering all these parameters, a larger sample size All rights reserved. No reuse allowed without permission.
(which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.
The copyright holder for this preprint this version posted May 2, 2023. ; 36 would be required. We expect such approaches to be feasible in the near future by performing metaanalyses across multiple (disease) population-based single-cell datasets that are uniformly processed in efforts such as those of the single-cell eQTLGen consortium (10).
In conclusion, this study highlights the importance of studying STEMI at cell-type-specific resolution, while taking genetic variation, disease status and phase of the disease into consideration. We expect that such an integrative approach as used here, will help to better grasp the molecular processes underlying STEMI and will be essential for the development of effective future therapies with reduced side-effects.

37
Processed (de-anonymized) scRNA-seq data, including a text file that links each cell barcode to its respective individual, has been deposited at the European Genome-Phenome Archive (EGA), which is hosted by the EBI and the CRG, under accession number EGAD00001010064.
All rights reserved. No reuse allowed without permission.
(which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.
The copyright holder for this preprint this version posted May 2, 2023. ; Figure S1. Overview of the study. Describing the cohort and the results per section All rights reserved. No reuse allowed without permission.

Supplementary figures
(which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.          (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

Supplementary tables
The copyright holder for this preprint this version posted May 2, 2023. (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.